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ACCELERATED CONVERGENCE FOR INCOMPRESSIBLE FLOW CALCULATIONS 


Geralyn M. Neely and Russell W. Claus 
National Aeronautics and Space Administration 
Lewis Research Center 
Cleveland, Ohio 44135 


Abstract 


Two improved a’gorithms which solve the 
steady-state Navier-Stokes equations, PISO and 
SIMPLER, are studied in this paper. Computations 
were carried out on progressively finer grids for 
the driven cavity and flow over a backward-facing 
step. The effects of relaxation factor, number of 
grid nodes and number of sweeps through the pres- 
sure equations are studied to evaluate the perfor- 
mance of the PISO and SIMPLER schemes. Results 
show that these improved schemes accelerate the 
convergence rate of the solution generally by a 
factor of two as compared to the SIMPLE method. 

Introduction 


Although currently available combustor models 
have already proven useful as analytical tools, 
these codes are not yet cost-effective alterna- 
tives to traditional design methods. One area in 
which combustor codes can be improved is in the 
execution time required for the solution algo- 
rithm. Typically, the method used to solve the 
steady-state Navier-Stokes equations is an itera- 
tive procedure, which, for stability reasons, must 
be underrelaxed. The widely-used SIMPLE (Semi- 
Implicit Method for Pressure-Linked Equations) 
scheme, developed by Patankar and Spalding in 
1972, is one such solution algorithm. * Although 
this method has been demonstrated to be quite 
effective, its convergence rate can be improved. 
This investigation focuses on two alternate 
approaches which a-celerate the solution to the 
steady-state Navier-Stokes equations. 

A short discussion of the solution algorithms 
studied, SIMPLE, SIMPLER (SIMPLE Revised!, 2 and 
PISO (Pressure Implicit Split Operator), ^ is in- 
cluded in this report. Results are presented for 
the driver, cavity and flow over a backwaro-facing 
step. The SIMPLE scheme is used as a basis for 
comparison to evaluate the performance of these 
algorithms. 


Symbol List 

C ,C ,C turbulence model constants 

u ci c 2 

k turbulence kinetic energy 

p pressure 

Re Reynold's number » 

uj fluctuating velocity about mean 

Ui mean velocity 

c turbulence energy dissipation rate 

f any of the independent variables 

o density 


turbulent viscosity 
turbulent Prandtl number 
Mathematical Model 


The TEACH code, developed at Imperial Col- 
lege, is the solution procedure used in this 
study.* This algorithm solves partial differen- 
tial equations which describe the fluid flow. 

Governing Equations 

The governing equations for two-dimensional, 
steady-state, incompressible, turbulent flow are 
shown below. The turbulence model employed is the 
two-equation k-c model developed by Jones and 
Launder. 5 

Continuity 



= 0 


Momentum 


u. U l!i_— ) 

1 3X . p 3X 3X. I p JX . 1 J ) 

— ^4 — v -7 v 


Convection Pressure 
gradient 


Diffusion 


Scalar transport 


U. 2SU.I ( JL- _ 77 \ 

i ax. p ax. o. ax. r 1 

1 v V 1 -t- - 1 U 


Convection 
Turbulence model 


Diffusion 


u i u j = 7 


3 U . 


Jx i 


P P \iX j 


u t = c P ^ 
t y c 


77— 1 u t 8? 

V “ ” p ^ 5x7 


j, UL..I j_piu_ — lUi 

i ax . p 3 x i \ a K 3 x i / 1 J 9x i 


U. = - -M — — ) - C 4 u-u • - C 
i a x i p ax^ l o c ax . j tj k i j 


ORIGINAL PAGE 19 
OF POOR QUALITY 


Solution Algorithms 

A staggered grid arrangement is set up to 
obtain the finite difference approximations.® 
Scalar variables such as pressure and turbulence 
kinetic energy are stored at mesh points, while 
velocities are stored at control volume faces. 
Integration is performed over the control volume 
surfaces. The partial differential equations are 
discretized using the hybrid (upwind/central) dif- 
ferencing scheme for the convective terms, and 
central differencing for the other terms. 

SIMPLE 

In the SIMPLE algorithm, a guessed pressure 
field is inserted into the discretized momentum 
equations to obtain a velocity field. The pres- 
sure field is corrected via an equation which is 
derived through a combination of the discretized 
continuity and momentum equations. The velocity 
field is then updated and is used to solve equa- 
tions for k, t and 9 . The corrected pressure 
field is treated as the guessed pressure field and 
the procedure is repeated until a converged solu- 
tion is obtained. 

The velocity correction equation used in the 
SIMPLE scheme for u at the point e is the 
fol lowing: 

y 

A u' = 'T A . U 1 . + (pn - Pc) A 

e e nb nb nb ' y ? H E' e 

Primes indicate corrections to old values. Tne 
underlined term, which represents the influence of 
corrected pressures on neighboring velocities, is 
neglected in the SIMPLE algorithm. To include 
this term would entail solving pressure correction 
equations for the whole flow field. On the other 
hand, the converged solution is unchanged by the 
exclusion of this term since this equation is 
merely used to correct the velocity. Doing this, 
however, results ir an overestimation of the pres- 
sure field, requiring relaxaton of the velocity 
equations to obtain a converged solution. Relaxa- 
tion factors less than 1.0 decelerate changes in 
the velocity from one iteration to the next. The 
lower the relaxation factor, the more slowly con- 
vergence is reached. The SIMPLE scheme necessi- 
tates the choice of very low relaxation factors to 
obtain convergence. 

SIMPLER 

The SIMPLER alqorithm improves upon the 
SIMPLE scheme by including the neglected terms in 
calculating the pressure field. The calculation 
sequence starts with a guessed velocity field. As 
in the SIMPLE scheme, an equation for the pressure 
field which combines the discretized momentum and 
continuity equations is solved. This pressure 
field is inserted into the discretized momentum 
equations to obtain a velocity field. A corrected 
pressure field is calculated and used to update 
the velocity field. This velocity field is 
treated as the guessed velocity field and the 
iteration procedure continues until convergence is 
reached. 

Because two additional equations are solved, 
each iteration step through the SIMPLER routine 
involves more computational time than an iteration 
step through the SIMPLE procedure. However, 


higher relaxation factors can be applied in the 
SIMPLER routine and therefore, convergence is much 
faster. 

PISO 

The PISO scheme also takes into account the 
terms neglected in the SIMPLE code, but in a dif- 
ferent manner. The PISO routine mimics the SIMPLE 
code up until the end of the first Iteration. At 
this point, the PISO scheme employs an equation 
containing the neglected terms which corrects the 
pressure. The velocity is then corrected and up- 
dated. Again, this procedure is repeated until 
the solution converges. For each iteration, the 
PISO scheme more closely satisfies continuity and 
momentum than the SIMPLE scheme, due to the addi- 
tional pressure and velocity corrections. In this 
manner, the PISO code allows for higher relaxation 
factors to accelerate convergence. 

Results 

The main thrust of this investigation is to 
analyze methods which speed the convergence of 
incompressible flow calculations. In order to 
accomplish this, several factors are considered. 
The most significant element in the optimization 
process is the relaxation factor. The choice of a 
too high relaxation factor may cause divergence; a 
too low relaxation factor results in a slow, 
costly calculation. Secondly, the optimum number 
of sweeps through the pressure equations per velo- 
city iteration is determined. This optimum varies 
according to the exactness with which each scheme 
solves the pressure-velocity equations. A third 
factor to consider is the grid size. A dense mesh 
may require much longer computational time than a 
coarse one. All of the above analyses are done 
for laminar flows. Finally, the "stiffness" im- 
posed by the turbulence model is studied. 

The two configurations chosen for this study 
are the square driven cavity and flow over a 
backward-facing step. The velocity vectors for a 
typical driven cavity calculation (Re = 1000) are 
illustrated (Fig. 1). The top wall of the cavity 
moves at a constant velocity, u, while the 3 other 
walls remain motionless. The grid for these cal- 
culations is evenly spaced. 

Figure 2 shows the geometry employed in the 
backward-facing step calculations. Computations 
for this geometry were performed at Reynolds num- 
bers of 1000 (laminar) and 30 000 (turbulent), 
based on the channel inlet height of 3 h (see 
Fig. 2). Aspect ratios (ax/ ay) for these calcula- 
tions range from 1.0 to 5.0. 

Effect of Relaxation Factor 

The importance of applying the optimum relax- 
ation factor for each scheme to obtain the fastest 
convergence time is illustrated (Fig. 3). Resid- 
ual error values are lxlO -5 for the laminar driven 
cavity and lxlO -4 for step flow calculations. 
Resiaual errors were chosen based on an overall 
error of 2 percent in the results. The minimum 
value on each curve represents the fastest conver- 
gence time and the corresponding optimum relaxa- 
tion factor. The vertical lines indicate the 
relaxation factors at which the calculations 
become unstable. For the 38 by 38 driven cavity 
calculations, (Fig. 3(a)), the SIMPLER and PISO 
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algorithms have higher optimum relaxation factors 
and therefore, lower convergence times. Note, 
however, that if a lower relaxation factor, for 
example 0.5, is chosen, both the SIMPLER and PISO 
schemes would be slower to converge than the 
SIMPLE routine. 

Similar results can be seen for the 62 by 22 
step flow calculations (Fig. 3(b)). For these 
computations, the PISO code has the lowest conver- 
gence time. Note also that with the PISO scheme, 
one can be quite a bit far from the optimum relax- 
ation factor with minimal effect. The opposite 
can be said of the SIMPLE algorithm, which has a 
very short range of relaxation factors which give 
reasonable results. Although the SIMPLER algo- 
rithm has a high optimum relaxation factor as com- 
pared to the SIMPLE code, its fastest convergence 
time is about equal to that of the SIMPLE code. 

Effect of Number of Pressure Sweeps 

Once the optimum relaxation factor for each 
scheme was established, the calculations were 
optimized with respect to the number of sweeps 
through the pressure equations. The optimum 
relaxation factor for each scheme is used in all 
further comparisons. The optimum number of pres- 
sure sweeps varies from 5 to 10, depending on the 
scheme and flow situation. Table I shows typical 
variations in convergence times for different num- 
bers of sweeps through the pressure equations for 
the PISO and SIMPLER algorithms. For the SIMPLE 
scheme, increasing the number of pressure sweeps 
beyond five resulted in an increased convergence 
time. Because of the neglected terms in the pres- 
sure correction equation, the pressure field at 
each iteration is quite inaccurate until the solu- 
tion nears convergence. Sweeping the inaccurate 
pressure field several times does little to 
enhance convergence for the SIMPLE algorithm. 

Since the PISO and SIMPLER schemes both include 
the neglected terms, the equations are more exact 
and the pressure field at each iteration is not 
nearly as inaccurate as that in the SIMPLE code. 

By increasing the number of sweeps through the 
pressure equations, the pressure field becomes 
more accurate and thus speeds convergence. 

Effect of Grid Size 


Computations were done for meshes ranging 
from 100 to 4700 nodes, with results shown in 
figure 4. All calculations were performed using 
the optimum relaxation factors and number of pres- 
sure sweeps determined previously. For the driven 
cavity (Fig. 4(a)) the PISO and SIMPLER schemes 
give nearly identical results, with the PISO 
scheme about 6 percent faster for the finest mesh. 
The relative efficiency, defined as 

((CPU time for SIMPLE) - (CPU for scheme)) 
(CPU for SIMPLE) 


remained at around 60 percent for both PISO and 
SIMPLER throughout the calculations. 

Laminar step flow results are shown in figure 
4(b). The SIMPLER scheme has the same convergence 
times as the SIMPLE code for grid sizes up to 1800 
nodes, and becomes faster for finer meshes, with 
an efficiency of about 40 percent at 4700 nodes. 


The PISO scheme efficiency improves from about 
20 percent for the coarse mesh to 40 percent for 
the 4700 node calculation. 

A convergence rate comparison (Fig. 5) gives 
insight into the directness or rapidity with which 
each scheme converges on a solution. These fig- 
ures also indicate that the results presented 
herein are not sensitive to the choice of conver- 
gence level. Again, optimum relaxation factors 
and number of pressure sweeps are used to obtain 
these results. For each iteration step, the maxi- 
mum residual error (made up of residuals of u 
momentum, v momentum and pressure) is established 
and graphed versus the execution time for that 
step. For the laminar driven cavity, the PISO and 
SIMPLER schemes have much steeper slopes and con- 
verge twice as fast as the SIMPLE algorithm for 
any given residual error value. In the step flow 
calculations, the PISO code again has a steeper 
slope and converges rapidly. The SIMPLE and 
SIMPLER codes performed similarly for this 
calculation. 

Turbulent Flow Calculations 


To investigate the "stiffness" imposed by the 
turbulence model, calculations were performed for 
the step flow geometry at a Reynold's number of 
30 000. "Stiffness" as it relates to turbulent 
flow calculations results from updatinq the turbu- 
lent viscosity at the end of each momentum 
pressure-correction iteration. In other words, 
the momentum equations are solved with the old 
values of turbulent viscosity and updated at the 
end of the iteration. This frequently slows con- 
vergence, i.e., imposes a "stiffness" on the cal- 
culation. These calculations (Fig. 6) were done 
for a residual error of 1x10 . The SIMPLER 

scheme gives the best performance, with the lowest 
convergence time at an optimum relaxation factor 
of 0 85 (Fig. 6(a)). The SIMPLER scheme also 
demonstrates a reasonable range of fast, stable 
operation for relaxation factors from about 0.80 
to 0.90. Although the PISO code improves upon the 
convergence time of the SIMPLE scheme by a factor 
of 2, it appears to have a much shorter range of 
stable operation than either of the other two 
schemes. 


For a series of increasingly finer meshes 
(Fig. 6(b)) the SIMPLER scheme speeds convergence 
by a factor of 3 with relative efficiencies rang- 
ing from 50 percent for the coarsest mesh to 60 
percent for the finest mesh. The PISO scheme is 

f iC f £»' f k ; r. M P ! F <- r* Kw -t f nm <•> ^ O tin. A 

I db Ltrr III OH LlIC 0 1 'If Ul ^UICI'IC u y u ■ ulvui t_ u<iu 

has a relative efficiency of 45 percent for all 
grid sizes. 


A convergence rate comparison for the turbu- 
lent flow case, again using the optimum relaxation 
factors for each scheme (Fig. 7), shows the PISO 
and SIMPLE schemes with steeper slopes and faster 
convergence times than the SIMPLER scheme for max- 
imum residual errors up to lx 10“ ^ . 


Conclusions and Recommendations 


Two improved solution algorithms, PISO and 
SIMPLER, have been tested in two flow situations, 
the driven cavity and flow over a backward-facing 
step. Effects of relaxation factor, number of 
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pressure sweeps and number of grid nodes on the 
convergence time of these schemes have beon ana- 
lyzed. These results show: 

1. The implementation of the PISO or SIMPLER 
scheme accelerates convergence generally by a fac- 
tor of 2 as compared to the SIMPLE algorithm. 

2. The efficiency relative to the SIMPLE 
scheme of the PISO and SIMPLER algorithms improves 
as more nodes are added to the flow field. 

3. The efficiency of the PISO and SIMPLER 
schemes is dependent upon the relaxation factor, 
flow situation and the number of sweeps through 
the pressure equations. 

Obviously, since these solution algorithms 
are tested here for two flow cases, these results 
cannot be taken as universal. However, they are 
useful in giving insight and making some general 
conclusions which should apply in similar flow 
situations. For a given class of problems, once 
the oot.mum relaxation factor and number of pres- 
sure sweeps has been determined, we recommend that 
these numbers be used for all further calculations 
in this class of problems. 

For general engineering use, it is recom- 
mended that a relaxation factor of 0.85 will give 
reasonable results for both the PISO and SIMPLER 
methods for laminar flow. The number of pressure 
sweeps to give fastest convergence is 10 for the 
SIMPLER method and 5 for the PISO method. For 
turbulent flow, a lower relaxation factor of 0.75 
for the PISO scheme and 0.85 for the SIMPLER 
scneme should yi'-e optimal results. Calculations 
for the SIMPLE method should be done with a relax- 
ation factor of 0.5 and 5 pressure sweeps for both 
laminar and turbulent flow. 


As a final note, it has been observed that 
the accuracy with which the pressure correction 
equations are solved affects the convergence time 
for the PISO and SIMPLER schemes. A faster tri- 
diagonal matrix algorithm solver may be imple- 
mented to make both of these schemes even more 
attractive. 
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TABLE I. DEPENDENCE OF CONVERGENCE TIME ON NUMBER 
OF PRESSURE SWEEPS 



Driven Cavity 38x38 
PISO 

No. pressure sweeps 

Conv. time, 1 


sec 

5 

74~2 

10 

104.1 

20 

136.0 

SIMPLER 

5 

94.9 

7 

83.6 

10 

78.5 

12 

85.2 


Laminar flow over a step 62x22 
PISO 

ReTaxat1cm~f actor I No. 


139.8 

139.1 

111.5 

151.4 


SIMPLER 
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Figure 1. - Velocity vectors for driven cavity calculations. 




RELAXATION FACTOR 
(b) Flow over a step. 

Figure 3. - Effect of relaxation factor on convergence time 
for laminar flow calculations. Reynolds number, 1000. 
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(b) Flow over a step. 

Figure 4. - Effect of grid size on convergence time for laminar 
flow calculations. 
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<b) Flow over a step. 

Figure 5. - Convergence rate comparison for laminar 
flow calculations. Reynolds number. 1000. 


CONVERGENCE TIME, sec CONVERGENCE TIME, sec 


ORfQfffi! 



(a) Convergence time versus rel< xation factor. 



(b) Convergence time versus number of grid points. 
Figure 6. - Turbulent flow calculations. 
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Figure 7. - Convergence rate comparison for turbulent 
flow calculations. Reynolds number, 30 00Q. 


